欧拉筛法筛积性函数

本文讨论的函数均为积性函数。

2.约数个数

方法一

网上大部分博客似乎都是这种方法,也很简单。

d[i]d[i]ii 的约数个数 , num[i]num[i]ii 的最小质因子的次数。

  1. ii 为质数

此时 d[i]=2d[i]=2num[i]=1num[i]=1

  1. prime[j]  iprime[j]~|~i

一个数的因数个数为:

i=1k(wi+1)   (wi为指数)\prod_{i=1}^k (w_i+1) ~~~ (\text{$w_i$为指数} )

prime[j]prime[j] 加入只会影响第一项,即让 w1+1w_1+1 变为 w1+2w_1+2

看一下定义,num[i]num[i] 就是 w1w_1 , 所以此时 d[iprime[j]]=d[i]/(num[i]+1)(num[i]+2)d[i*prime[j]]=d[i]/(num[i]+1)*(num[i]+2)

因为 prime[j]prime[j]iprime[j]i*prime[j] 的最小质因子,所以 num[iprime[j]]=num[i]+1num[i*prime[j]]=num[i]+1

  1. prime[j]  iprime[j]~ \not|~i

由积性函数的性质得:d[iprime[j]]=d[i]d[prime[j]]d[i*prime[j]]=d[i]*d[prime[j]]

因为 prime[j]prime[j]iprime[j]i*prime[j] 的最小质因子,所以 num[iprime[j]]=1num[i*prime[j]]=1

void sieve( ) {
    d[ 1 ] = 1;
	for( int i = 2 ; i <= MAXN ; i ++ ) {
        if( !vis[ i ] ) {
            prime[ ++ k ] = i;
            d[ i ] = 2 , num[ i ] = 1;
        }
        for( int j = 1 ; j <= k && 1ll * i * prime[ j ] <= MAXN ; j ++ ) {
            vis[ i * prime[ j ] ] = 1;
            if( i % prime[ j ] == 0 ) {
                d[ i * prime[ j ] ] = d[ i ] / ( num[ i ] + 1 ) * ( num[ i ] + 2 );
                num[ i * prime[ j ] ] = num[ i ] + 1;
                break;
            }
            d[ i * prime[ j ] ] = d[ i ] * d[ prime[ j ] ];
            num[ i * prime[ j ] ] = 1;
        }
    }
}

方法二.

上面的方法需要多开一个数组,太浪费空间了,介绍一种更简单的方法。

只讨论 prime[j]  iprime[j]~ |~ i 的情况,其余同上。

由约数个数公式得:

d[i]=j=1k(wj+1)=(w1+1)×j=2k(wj+1)d[i]=\prod_{j=1}^k(w_j+1)=(w_1+1) \times \prod_{j=2}^k(w_j+1)

因为后面的根 prime[j](p1)prime[j](p_1) 无关,为了方便记为 tt

又有:

d[i×p1]=(w1+1+1)×t=d[i]+td[i \times p_1]=(w_1+1+1) \times t = d[i]+t

d[i/p1]=(w1+11)×t=d[i]td[i/p_1]=(w_1+1-1) \times t = d[i]-t

两式相加得:

d[i×p1]+d[i/p1]=2×d[i]d[i \times p_1]+d[i / p_1]=2 \times d[i]

移相得:

d[i×p1]=2×d[i]d[i/p1]d[i \times p_1]=2 \times d[i]-d[i/p_1]

d[i×prime[j]]=2×d[i]d[i/prime[j]]d[i \times prime[j]]=2 \times d[i]-d[i/prime[j]]

void sieve( ) {
    d[ 1 ] = 1;
    for( int i = 2 ; i <= MAXN ; i ++ ) {
        if( !vis[ i ] ) {
            prime[ ++ k ] = i;
            d[ i ] = 2;
        }
        for( int j = 1 ; j <= k && 1ll * i * prime[ j ] <= MAXN ; j ++ ) {
            vis[ i * prime[ j ] ] = 1;
            if( i % prime[ j ] == 0 ) {
                d[ i * prime[ j ] ] = 2 * d[ i ] - d[ i / prime[ j ] ];
                break;
            }
            d[ i * prime[ j ] ] = d[ i ] * d[ prime[ j ] ];
        }
    }
}

3.约数和

f[i]f[i]ii 的约数和。

  1. ii 为质数

此时 f[i]=i+1f[i]=i+1

  1. prime[j]  iprime[j]~|~i

由约数和公式有:

f[i]=d=1kj=0wdpdj=j=0w1p1j×d=2kj=0wdpdj=d=0w1p1d×tf[i]=\prod_{d=1}^k\sum_{j=0}^{w_d}{p_d}^j=\sum_{j=0}^{w_1}{p_1}^j \times \prod_{d=2}^k\sum_{j=0}^{w_d}{p_d}^j=\sum_{d=0}^{w_1}{p_1}^d \times t

其中 pdp_d 为第 dd 个质数 , wdw_d 为第 dd 个质数的指数。

prime[j]prime[j] 的加入只会影响含 p1p_1(最小质因子) 的式子,所以后面的根本不用考虑,为了方便记为 tt。因为 prime[j]prime[j] 就是 p1p_1 , 下面简写。

考虑 f(i×p1)f(i \times p_1)f(i/p1)f(i/p_1)

f(i×p1)=d=0w1+1p1d×t=f(i)+p1w1+1×t            (1)f(i \times p_1)=\sum_{d=0}^{w_1+1}{p_1}^d \times t= f(i)+p_1^{w1+1} \times t ~~~~~~~~~~~~ (1)

f(i/p1)=d=0w11p1d×t=f(i)p1w1×tf(i / p_1)=\sum_{d=0}^{w_1-1}{p_1}^d \times t= f(i)-p_1^{w1} \times t

p1×f(i/p1)=d=0w11p1d×t=p1×f(i)p1w1+1×t    (2)\Rightarrow p_1 \times f(i / p_1)=\sum_{d=0}^{w_1-1}{p_1}^d \times t= p_1 \times f(i)-p_1^{w1+1} \times t ~~~~ (2)

(1)+(2)(1)+(2) 得:

f(i×p1)+p1×f(i/p1)=(p1+1)×f(i)f(i \times p_1)+ p_1 \times f(i / p_1)=(p_1+1) \times f(i)

f(i×p1)=(p1+1)×f(i)p1×f(i/p1)f(i \times p_1)=(p_1+1) \times f(i)-p_1 \times f(i/p_1)

还原得:

f(i×prime[j])=(prime[j]+1)×f(i)prime[j]×f(i/prime[j])f(i \times prime[j])=(prime[j]+1) \times f(i) - prime[j] \times f(i/prime[j])

  1. prime[j]  iprime[j]~ \not|~i

由积性函数的性质得:f[iprime[j]]=f[i]f[prime[j]]f[i*prime[j]]=f[i]*f[prime[j]]

void sieve( ) {
    f[ 1 ] = 1;
    for( int i = 2 ; i <= MAXN ; i ++ ) {
        if( !vis[ i ] ) {
            prime[ ++ k ] = i;
            f[ i ] = i + 1;
        }
        for( int j = 1 ; j <= k && 1ll * i * prime[ j ] <= MAXN ; j ++ ) {
            vis[ i * prime[ j ] ] = 1;
            if( i % prime[ j ] == 0 ) {
                f[ i * prime[ j ] ] = ( prime[ j ] + 1 ) * f[ i ] - prime[ j ] * f[ i / prime[ j ] ];
                break;
            }
            f[ i * prime[ j ] ] = f[ i ] * f[ prime[ j ] ];
        }
    }
}